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It has proved difficult to extend tlie density matrix renormalization group technique to large two- 
dimensional systems. In this Communication I present a novel approach where the calculation is 
done directly in two dimensions. This makes it possible to use an infinite system method, and for 
the first time the fixed point in two dimensions is studied. By analyzing several related blocking 
schemes I find that there exists an algorithm for which the local energy decreases monotonically as 
the system size increases, thereby showing the potential feasibility of this method. 



There is a great need for improved numerical tech- 
niques that are able to treat two-dimensional electronic 
systems. Fundamental questions such as the existence 
of superconductivity in the 2D Hubbard and t-J mod- 
els have not been resolved with present-day analytical or 
numerical techniques. All the commonly used numerical 
techniques suffer from shortcomings: exact diagonaliza- 
tion is limited to small lattice sizes duei-to the exponential 
growth of states with the system size.tl Quantum Monte 
Carlo calculations are unable to reach large fermion sy, 
tems at low temperatures due to the "sign problem".! 
DMRG calculations have be^ very successful at treat- 
ing one-dimensional systems ,0 but accurate results are 
difficult to obtain for large two-dimensional systems. 

It seems most likely that major improvements towards 
a reliable 2D technique will be made within the DMRG 
context. Years of effort have resulted in iip progress to- 
wards solving the fermion "sign problem" .H It also seems 
unlikely that the size of available computer memory will 
increase fast enough to facilitate exact diagonalizatio; 
calculations for large systems. Recent DMRG studiesjS'l 
on the other hand, have reached the largest 2D systems 
reported to date. In this Communication I approach the 
2D DMRG calculation from a new angle, which hopefully 
may encourage further research in this direction. 

First the basis of standard ID DMRG and previous 
2D DMRG calculations will be reviewed. The source of 
the difficulties with previous 2D calculations is discussed. 
Thereafter several new 2D blocking schemes are proposed 
and tested, keeping only a small number of states per 
block. Finally a promising algorithm is discussed in more 
detail. All calculations are performed on the 2D antifer- 
romagnetic Heisenberg model. The ground state param- 
eters of this fundamental auantum-mechanical model are 
known to a high accuracy)] which makes it an ideal test- 
ing ground for new numerical techniques. 

The central idea in a DMRG calculation is to itera- 
tively increase the system size, but to avoid an exponen- 
tially increasing number of states by keeping only a fixed 
number of the "most important" states at each iteration. 
In early numerical renormalization calculations the low- 
est eigenstates of the system were chosen to be the. "im- 
portant states" , but the results were discouraging.! The 
major breakthrough came with White's insightB to use 



the density matrix to determine which states to keep. In 
the superblock method a number of blocks are combined 
together to form a superblock. The superblock is divided 
into two parts, the "system block" and the "environment 
block" . At each iteration the superblock is diagonalized 
and the density matrix is formed for the system block. 
The density matrix is diagonalized and the importance of 
each eigenstate is given by its eigenvalue. The states with 
the largest eigenvalues are kept and the rest discarded. 
In the next step of the iteration this system block, with 
a reduced number of states, will be used in forming the 
new superblock. Thus the number of sites increases with 
each iteration, while the number of states kept remains 
constant. 

Using this basic formula many different DMRG algo- 
rithms can be defined. An algorithm is characterized by 
how the superblock is constructed, and by the manner 
in which the blocks are enlarged. The most commonly 
used method was proposed in White's original work. The 
superblock consists of four blocks, with the two central 
blocks consisting of single sites, and the two end blocks 
being reflections of each other, see Fig. 0. The system 
block is taken to be half the superblock, that is, one end 
block and an adjacent single site, here called a building 
block. Using the density matrix a fixed number of states 
are kept for the system block, which in the next itera- 
tion will be recombined with a building block to form a 
new system block. In this manner the system blocks are 
enlarged by the size of the building block (here one site) 
at each iteration. A variation of this method is to form 
the superblock out of three blocks, two end blocks and a 
single site in the middle. The system block is chosen as 
above, but the environment block consists of a single end 
block. 

If the above procedure is iterated repeatedly one can 
reach arbitrarily large systems. This is called the infi- 
nite system method. In the thermodynamic limit the 
energy does approach a fixed point, and as the number 
of states kept is increased, the fixed point approaches the 
bulk ground state energy for the model. Usually the in- 
finite system method has been used to measure various 
quantities in the middle of very large systems. In this 
manner results-iwith an accuracy of up to 13 digits have 
been reportedJj 




systems to date, up to sizes 11 by 16fl An alternative 
approach is to add a row of sites at each iteration. In 
this manner strips with a width of jjp to 6 sites and a 
length of 30 sites have been studied^ 



FIG. 1. Standard ID blocking algorithm. The superblock 
consists of the system block (heavy dashed line) and the envi- 
ronment block (dashed line). The building block is indicated 
by a dotted line. 

The rate of convergence does depend on how the su- 
perblock is constructed. If one tries to add the building 
block at the outer boundary of the superblock, then the 
local energy will increase as one increases the system size, 
indicating that a good basis has not been chosen. This 
can easily be understood since with open periodic bound- 
ary conditions the wave function has to vanish at the 
boundary. If one adds a site to the boundary it is clear 
that some artifact will remain in the wave function as the 
system size is increased. This somewhat trivial example 
shows that one cannot construct an arbitrary blocking 
scheme and expect that a DMRG calculation will yield a 
fast convergence. This will become more evident for 2D 
systems. 

In addition to the infinite system method discussed 
above, a finite system method, also introduced in White's 
pioneering work, is commonly used. This algorithm is 
used to calculate properties of finite size systems to high 
accuracy. Initially the infinite system method is used to 
reach a system of desired length L. At each iteration the 
system block is saved, so that when a system of size L 
is reached, system blocks of sizes 1 to L/2 — 1 are saved. 
Once the desired total system size has been reached the 
superblock size is fixed. Next the saved system blocks 
are used as environmental blocks while the system block 
size is increased until it has reached the maximum length 
L — 3. Now system blocks of sizes 1 to L — 3 are saved 
and these blocks can be used as environment blocks for 
a consequent sweep through the lattice. In this manner 
the basis kept in the system blocks can be iteratively 
improved until convergence is reached. 

After this brief review of ID DMRG calculations previ- 
ous 2D calculations will be considered next. Most previ- 
ous 2D calculations involve mapping thCrSD lattice onto a 
ID system with long-range interactions ,EHl3~t3 see Fig. 0. 
Thereafter the above ID finite system algorithm is used. 
Notice that using this mapping it is not possible to use 
an infinite system algorithm, since one determines the 
size of the final lattice when doing the mapping. Also, 
the blocks break the symmetry of the lattice and it is 
generally not possible to use a reflection of the system 
block as the environment block. Therefore one has to 
use some different trick to form the environment block for 
the initial sweep. The two simplest options are to either 
leave the environment block empty, or set all long-range 
interactions to zero in the initial sweep. This method 
has, however, been able to treat the largest 2D fermion 




FIG. 2. Mapping of the 2D lattice to a ID lattice with 
long-range interactions, the solid line indicates the Id system, 
while the dashed lines show the long-range interactions. 

Why have larger systems not been studied ? Liang 
and PangiJ found that for a 2D gas of free electrons, the 
number of states needed to maintain a certain accuracy 
grows exponentially with the linear system size. This 
convergence was also confirmed for an algorithm were a 
row of sites was added at each step.El Although no proof 
has been given, this statement is often referred to as most 
probably valid for any 2D DMRG calculation. This state- 
ment was, however, made for small finite size systems and 
it is not clear that it will apply to possible infinite sys- 
tem methods. In an infinite system DMRG calculation a 
fixed number of states are kept as the size of the system 
is increased. According to the above statement accuracy 
should be lost in the process. For a system with open 
boundary conditions the local energy decreases as the 
system size is increased. Furthermore, due to the varia- 
tional character of the technique,tj the DMRG energy is 
an upper bound on the energy of the system. Therefore 
accuracy would certainly be lost if the DMRG energy in- 
creased as the system size is increased, in agreement with 
the above statement. But if the energy decreased as the 
system size is increased, then the bound on the system 
energy is continuously improved, and in the limit of the 
fixed point the relative accuracy will approach a constant 
although only a fixed number of states is kept. 

It was therefore the goal of this study to investigate 
whether there exist 2D blocking algorithms for which 
the energy decreases monotonically as the system size 
is increased. The algorithm should retain more of the 
symmetry of the lattice so that it can be used in the in- 
finite system mode, and the fixed point studied directly. 
The reason for the above mapping of the 2D system to 
a ID system is that it is not trivial to construct such an 
algorithm. Since the superblock of most symmetric 2D 
algorithms is bound to consist of more blocks than in the 
ID case computer memory limitations will also be more 
severe. 

In order to build up a two-dimensional lattice in a more 
symmetric fashion it seems likely that one has to use 
building blocks consisting of several sites. The simplest 



idea is probably io use a row of sites as building block 
at each iteration.Q This method has, however, two short- 
comings; it only grows the lattice in one spatial dimen- 
sion, and the number of states in the added row increases 
exponentially with the width of the lattice. 

In a first attempt to overcome these problems I di- 
vided the square 2D lattice up into three blocks, consist- 
ing of the diagonal, a triangular block below the diagonal, 
and the reflection of this block above the diagonal, see 
Fig. H The lower triangular block is used as the sys- 
tem block, the diagonal as the building block and the 
reflection of the lower triangular block is used as the en- 
vironment block. At each iteration the whole diagonal is 
thus added to the lower triangular block. In this manner 
one of the problems with adding just a row of sites to 
the system block is overcome; the lattice grows in both 
spatial dimensions. Furthermore, the procedure retains a 
high degree of symmetry and could, in principle, be used 
as an inflnite system method. The problem is, of course, 
that the number of states needed to describe the exact di- 
agonal block still increases exponentially with the linear 
system size. The method was, however, implemented. 





FIG. 3. A blocking algorithm for the 2D lattice where the 
system block (heavy dashed line) is formed by combining the 
lower triangular part of the lattice with the diagonal block ( 
dotted line) . The upper triangular part of the lattice is used 
as environment block (light dashed line). 

When adding an exact diagonal to the system the en- 
ergy per site decreased as the system size was increased, 
until the computer ran out of memory. Having passed 
this simple test the next problem to be addressed was 
the exponential increase of states in the diagonal. The 
natural way to avoid the exponential growth is to select 
only the most important states in the diagonal block by 
diagonalizing the density matrix for the diagonal. This 
was done, and at each iteration a single site was added 
in one corner. The local energy did, however, start to 
increase as the system size was increased. As in the ID 
case, the reason seemed to be that a site was added at the 
boundary of the system. If periodic boundary conditions 
are used this may be a possible blocking formula, but it 
does not work with open boundary conditions. 

A new method was therefore tried where, in analogy 
with the ID method, the diagonal was divided into two 
blocks with the additional site added in the middle. The 
local energy did, however, still increase as a function of 
lattice size. A potential problem seemed to be that when 
working with square lattices, one is forced to construct 



the wave function for a lattice with an even number of 
sites from the wave function for a lattice with an odd 
number of sites. Lattices with odd and even numbers 
of sites do, however, have quite different wave functions. 
This issue can be avoided if one studies lattices tilted 
by 45 degrees. First an attempt was made to use lat- 
tices tilted by 45 degrees containing an even number of 
sites, see Fig. ^ With this geometry the standard ID 
DMRG technique can be used for the diagonal, with the 
exception that the superblock also contains the triangular 
blocks. The site energy did still not decrease monoton- 
ically as the system size was increased. Therefore tilted 
lattices with an odd number of sites were investigated, 
see Fig. ||, and it was found that the energy decreased 
monotonically as the lattice size was increased. This was 
the only blocking scheme found in this study for which 
this was the case. The fact that there exists such an 
algorithm is certainly encouraging, and not self-evident, 
as pointed out above. Since this was the most promising 
algorithm found in this study the results will be analyzed 
in more detail next. 
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FIG. 4. A blocking algorithm for a tilted 2D lattice with 
an even number of sites. The density matrix is formed both 
for the diagonal block and the upper triangular block (both 
indicated by heavy dashed lines). 





FIG. 5. A blocking algorithm for a tilted 2D lattice with 
an odd number of sites. The density matrix is formed both 
for the diagonal block and the upper triangular block (both 
indicated by heavy dashed lines). 

In Fig. |6| the site energy is shown as a function of the 
number of iterations. Density matrices are formed both 
for the diagonal and the triangular blocks, and the num- 
ber of states kept in these blocks are denoted rud and 
rrit respectively. The ground staije energy for the 2D 
Heisenberg model is -0.669437(5).Q Keeping four states 



in the diagonal and four blocks in the triangular block 
the energy levels out around -0.49. Increasing m^ to 16 
dramatically improves the energy to -0.57. It seems de- 
sirable to keep a higher number of states in the diagonal 
than in the triangular block. Next the number of states 
in the triangular block was increased to eight. Then only 
16 states could be kept in the triangular block, and for 
the number of iterations that could be done the results 
were slightly better than the results obtained when keep- 
ing four states in the triangular block and 16 states in 
the diagonal block. 
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FIG. 6. The energy as a function of number of iterations 
for the 2D Heisenberg model. The number of states kept in 
the diagonal block is denoted irid, while the number of states 
in the triangular block is given by rrit. The exact result is 
obtained from Monte Carlo simulations (Ref. 7). 

In order for a DMRG method to be useful one has to be 
able to keep enough states per block to reach convergence 
in the quantity studied. Since a 2D blocking algorithm 
of the kind described above contains more blocks than 
the traditional ID blocking method this may prove difh- 
cult. The superblock for the above algorithm will contain 
2rri^frif states, and the density matrix for the triangular 
block will contain 2m^mt states. The programs used in 
this investigation are, however, far from optimized. Com- 
puter memory limited the present study. By using good 
quantum numbers, like the z-component of the spin for 
the Heisenberg model, all matrices become block diag- 
onal. In this study complete matrices were stored, and 
one should be able to significantly increase the number of 
states kept if only the non-zero matrix blocks are stored. 

The important issue is thus to study how the bulk 
ground state energy is approached as one further in- 
creases the number of states in the blocks. The reason for 
the great success of the ID DMRG method is a very fast 
convergence. If the 2D calculation shows exponential or 



power law convergence one may be able to keep enough 
states to reach accurate results, but if the convergence is 
slower this may prove difficult. 

It is also possible to make a finite system algorithm 
based on the above blocking procedure. In such a method 
one could use the infinite system method for the initial 
steps, saving both the triangular and the diagonal blocks. 
Further sweeps could use these saved blocks as environ- 
ment blocks and improve the basis kept in the system 
blocks, in a manner analogous to the ID method. 

The algorithm presented in this Communication bears 
some resemblance to a "four-block method" proposed by 
Bursill.Efl In both methods the building block, which de- 
termines the growth of the system block, does not consist 
of exact sites, as in the original method, but of sites with 
a reduced number of states. 

To conclude, using a new approach I have explored the 
first fully two-dimensional infinite system DMRG calcu- 
lation. The fixed point in two-dimensions could be ex- 
plicitly studied and it was shown that there exists an 
algorithm for which the local energy for the Heisenberg 
model decreases monotonically as the system size is in- 
creased. Previous results indicated that it is necessary to 
keep a number of states that grows exponentially with the 
linear systein size to maintain a certain accuracy. This 
does not appear to be the case with the infinite system 
algorithm as the fixed point is approached. The method 
preserves a high degree of the symmetry of the lattice and 
could be used as a starting point for a finite system al- 
gorithm. Further studies are necessary to verify whether 
the method presented here, or other similar algorithms, 
exhibit convergence that is fast enough to calculate prop- 
erties of large two-dimensional electronic systems. 
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